Orthogonal versus covariant Lyapunov vectors for 
rough hard disk systems 



H. Bosettfj and H.A. Posch 

Computational Physics Group, Faculty of Physics, University of Vienna, 
Boltzmanngasse 5, 1090 Wien, Austria 

E-mail: hadrien . bosett i@univie . ac . at 

E-mail: harald . pos chOunivie . ac . at 

Abstract. The Oscledec splitting of the tangent space into covariant subspaces 
for a hyperbolic dynamical system is numerically accessible by computing the full 
set of covariant Lyapunov vectors. In this paper, the covariant Lyapunov vectors, 
the orthogonal Gram-Schmidt vectors, and the corresponding local (time-dependent) 
Lyapunov exponents, are analyzed for a planar system of rough hard disks (RHDS). 
These results are compared to respective results for a smooth-hard-disk system 
(SHDS). We find that the rotation of the disks deeply affects the Oseledec splitting 
and the structure of the tangent space. For both the smooth and rough hard disks, 
the stable, unstable and central manifolds are transverse to each other, although the 
minimal angle between the unstable and stable manifolds of the RHDS typically is 
very small. Both systems are hyperbolic. However, the central manifold is precisely 
orthogonal to the rest of the tangent space only for the smooth-particle case and not 
for the rough disks. We also demonstrate that the rotations destroy the Hamiltonian 
character for the rough-hard-disk system. 
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1. Introduction 

The phase-space trajectory of particle systems generally is sensitive to tiny 
(infinitesimal) perturbations. The ensuing exponential instability is revealed by 
following the time evolution of a complete set of perturbation vectors, {5Ti}® =1 , spanning 
the tangent space, and by computing the whole set of exponential rate constants, 
the Lyapunov exponents. Here, D denotes the phase-space dimension. Two kinds of 
perturbation vectors are usually employed: the orthonormal Gram-Schmidt vectors, 
{<7j}° =1 , and the generally non-orthogonal covariant Lyapunov vectors, {"Uj}? =1 . The 
former were introduced in 1980 by Benettin et al. [U [2] and by Shimada and Nagashima 
[3] for the computation of the spectrum of Lyapunov exponents. They are made 
orthonormal by a periodic applications of a Gram-Schmidt procedure to counteract 
the tendency of all perturbation vectors to rotate towards the most unstable direction 
in phase space. This procedure is equivalent to a QR decomposition at periodic time 
intervals jU [5]. For a recent review we refer to Ref. [6]. The covariant vectors, on 
the other hand, are generally non-orthogonal and were made numerically accessible in 
2007 by Ginelli et al. |7j and by Wolfe and Samelson [8]. They co-move (co-reorient, 
in particular) with the flow in tangent space and define a hierarchy of directions at 
any phase-space point, which are individually associated with the Lyapunov exponents 
(as long as the exponents are non-degenerate). The algorithm for their construction is 
more involved and still requires the computation of the Gram-Schmidt vectors in the 
forward direction of time, followed by an iteration in the backward direction to obtain 
the covariant vectors (see [9j [10]). Both vector sets are intimately related to various 
subspaces of the tangent space introduced by Oseledec in 1968 [TT| . The Gram-Schmidt 
vectors are related to the eigenvectors of the Oseledec matrix j5j [12]. The covariant 
vectors correspond to the Oseledec splitting [131 EH M an d> thus, provide a very natural 
spanning set of the tangent space at any phase point. 

The time-averaged growth rate of a covariant vector gives rise to a (global) 
Lyapunov exponent, and the set of exponents, {Ai}° =1 , ordered according to size, 
Aj > A i+ i, is referred to as the Lyapunov spectrum [IT] . The sub-exponential growth 
rates give vanishing exponents, for which the associated covariant vectors span the 
central manifold. The covariant vectors with positive (negative) Lyapunov exponents 
span the unstable (stable) manifold. For conservative systems the sum of all Lyapunov 
exponents vanishes. This sum is negative for dissipative systems, which implies a 
collapse of the phase-space distribution onto a multifractal strange attractor |T5] . At 
any phase-space point local Lyapunov exponents may be defined. They vary with time 
when the system evolves. The global exponents are time averages of the local Lyapunov 
exponents [TH [17] . 

In this paper we study the Lyapunov instability for rough-hard-disk systems 
(RHDS) [18]. By adding the notion of roughness and angular momentum to each 
particle, an exchange of energy between translational and rotational degrees of freedom 
is introduced. The three-dimensional version of this model - arguably - may be 
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considered as the simplest model of a molecular fluid. Our results are compared to those 
for smooth-hard-disk systems (SHDS), which have been the subject of many previous 
studies pll EH EH E2J E3] . 

The paper is organized as follows. The basic concepts of the dynamics of the 
tangent-space bundle are summarized in Sec. EJ The smooth-hard-disk system and, 
more thoroughly, the rough-hard-disk system are introduced in Sec. [31 In Sec. H] we 
compare the Gram-Schmidt vectors to the covariant Lyapunov vectors (which give rise 
to the same global Lyapunov exponents) for the RHDS. Our aim is to characterize the 
transition from non-rotating (smooth) to rotating (rough) particles. The null subspace 
for the RHDS is the topic of Sec. El where we show that it is not strictly orthogonal 
to the unstable and stable subspaces. Using the covariant vector set, we show that the 
dynamics is converging to a sub-manifold. In Sec. [6] we examine the loss of isolation for 
the local covariant exponents, which implies that a finite-time average of local exponents 
may lead to a different order in the Lyapunov spectrum as compared to that of globally- 
averaged exponents. Furthermore, the transversality between the stable, the unstable, 
and the central manifolds is investigated for the SHDS and the RHDS. Sec. [7] is devoted 
to a discussion of the qualitative differences between the two models and of the break 
down of the Hamiltonian character for the RHDS due to the coupling between the 
translational and rotational degrees of freedom. We conclude in Sec. [8] with a short 
summary. 

2. Tangent-bundle dynamics 

First we recall some basic notions in statistical mechanics and dynamical systems theory 
(see e.g. [21]): (i) The concept of a phase space X is tantamount to consider the full 
ensemble of accessible system states, denoted by the state vector T. If this system is 
ergodic, which we assume, every infinitesimal volume element of the phase space in 
accord with the conserved quantities will be visited in due time. Each phase point T 
on the hypersurface, defined by the exclusion of the subspace violating the conservation 
laws, becomes equally likely; (ii) The dynamics in the phase space is characterized by 
the motion of a phase point in that space with a fixed frame, which is not the case for 
the tangent space TXp. TXp moves with the state vector in phase space and is tangent 
to the hypersurface determined by the conserved quantities. It contains the complete 
set of perturbation vectors which point in all accessible directions, where the state 

r is allowed to evolve (observing the conservation rules); (Hi) In order to examine also 
the dynamics of the perturbation vectors ST, we have to consider the concept of the 
tangent bundle TX, which consists of the union of the tangent spaces to X at various 
state vectors T, [J reX TXr. A point of TX is a vector ST, tangent to X at T(t). 

According to the multiplicative ergodic theorem of Oseledec for unstable systems 
[TT] . the dynamics of the tangent bundle induces some hierarchically dissociated sub- 
bundles. At some phase point T, the subbundles correspond to a set of dissociated 
subspaces, which is referred to as Oseledec splitting. The exponential growth (or decay) 
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of the subbundles gives rise to the Lyapunov spectrum, whose hierarchical structure 
is reflected by the ordering of the exponents by size. Ginelli et al. [7] introduced an 
algorithm which permits to generate numerically the Oseledec splitting via the covariant 
vectors {«,:}. This method is based upon a theoretical statement by Ruelle [13] and 
Ershov and Potapov [12] . 

Vi = (g{ +) © ■ ■ ■ © g[ +) ) n (g^ © • ■ ■ © g^ , (l) 

where Vj is a one-dimensional subspace of TXp spanned by the covariant vector Vi, 
and denotes the set of orthonormal Gram-Schmidt vectors. As usual, D is the 
phase space dimension. All of these vectors are associated with the same state vector 
r. The superscript "+" refers to the forward direction of time (i.e. the starting vector 
set {g\ } was initially placed in the far past), whereas the superscript "— " refers to 
the backward direction of time (i.e. the starting vector set {g\ } was initially placed 
in the far future). 

Let T(t) denote the state of the system at time t. The flow in phase space is denoted 
by 0* : X — > X and the corresponding flow in tangent space by Dcp 1 : TX r ( ) — > TXr(t). 
Then the covariant vectors {v{\ obey 

iV|r(o)<r(o)) = «<(r(t)) , 

for all i G {1,...,D}, where T(£) = 0'(r(O)). The corresponding (global) Lyapunov 
exponents are time averages, 

lim 1 log || D^\ m Vi(T(0)) || = ±Xi , 

and due to the dynamical hierarchy, they obey: Ai > ■ • • > A D . When we demonstrate 
covariant vectors in our numerical work below, we shall always consider them as 
normalized (without explicit reference to the normalization procedure). 



3. Systems of rough hard disks (RHDS) 

Rough disks temporarily store energy in internal degrees of freedom (rotation). They 
suffer elastic hard collisions and move along straight lines inbetween. The systems 
consist of N particles at a density p in a box of dimension (L x , L y ) with periodic 
boundaries in the x- and y-directions. The state vector is given by 

T = {<fn, Pn,^n}n=l > ( 2 ) 

where q n , p n ,u) n denote the respective position, momentum and angular velocity of 
particle n. The orientation angles of the disks are not required for the equations of 
motion and are not included in the list of independent variables. This means that the 
Lagrangian of the system does not depend on the orientation of the particles. We shall 
come back to this point in Sec. El An arbitrary perturbation vector is given by 

ST = {5q n , 5p n , 5u n }% =1 . (3) 
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The rough hard-disk system (RHDS) is an extension of the common smooth 
hard-disk system (SHDS), where the particles rotate and exchange translational and 
rotational energy at a collision. This implies a coupling between translational and 
rotational degrees of freedom. The "rough-hard-collision" transformation rule is 
established from the conservation laws for energy, linear momentum, and angular 
momentum. A first version of the RHDS was already introduced in the late 19 
century by Bryan [25]. Explicit motion equations were given by Chapman and Cowling 
[26] in 1939. Initially, roughness was introduced as the maximum possible roughness 
in the model, where a collision reverses the relative surface velocity at the point of 
contact of two colliding particles. Such a model was studied by O'Dell and Berne 
[27]. Subsequently, Berne and Pangali also investigated models with partial roughness 
[28] [29] . Here we use the maximum-roughness model of Chapman and Cowling [26J. 

The energy of the RHDS is purely kinetic, 

1 - 



m 



where I is the moment of inertia, which is taken to be identical for all disks. 

We first recall the collision rules for two rough hard spheres k and i in three 
dimensions [26], [18]. Let q = q k — qe and Sq = 5qk — 8qe denote the respective relative 
positions in the phase space and in the tangent space, and analogously, p = Pk — Pe 
and dp = 5pk — o~pe the relative momenta. For convenience, we define the following 
quantities: Q = Uk + tie, 6Q = 5uk + Sue, the unit vector n — q/a, the dimensionless 
moment of inertia k = 4/ /ma 2 , and the velocity of the impact point on the surface of 
the scatterer, 

9 = - + -(nxn) . 
m 2 

During a collision, only the components associated with the colliding particles k and I 
of the vector Y change [IB]. 

9* = <ffc , 
Qt = & , 

fk =Pk-m<yg- (5{n-p)n , 

Pi = Pi + mi g + 0{n ■ p)n , 

-* 2/3 _ . 

uj k = uj k H (n x g) , 

a 

2/3 

3 1 = <2 t -\ (nx g) , (4) 

a 

where the prime, here and below, refers to the state immediately after the collision. The 
constant parameters 7 = + ac) and (3 = 1/(1 + k) control the coupling between 
translational and rotational degrees of freedom. If the moment of inertia vanishes (k — > 
and f3 — > 1) the translational and rotational degrees of freedom decouple. 
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Using the general equation, Eq. (14) of Ref. 
the collision map for the perturbation vectors 5F 



the linearization of Eq. (J3J) yields 



5q e - 



70 H (n • p)n 

m 

7^ H (n ■ p)n 



m 



5r r 



5r r 



ft 



5pk — m^5g c [(5q c • p)n+ cr(n ■ 5p)n + (n ■ p)5q c 



Sty 



5(2! 



Sp e + rwySg, 
2/3 



a 

ft 



52 k + 

but + ~r 
a 1 



— [(8q c ■ p)n + a(n ■ 5p)n + (n ■ p)5q c 
a 

Sq c x g + a(n x 
8q c x 9 + o-(n x 5<f c 



(5) 



where 5t c = —(5q ■ n)/(p/m ■ ft) is the infinitesimal time shift between the collision of 
the reference and the perturbed trajectories, and Sq c = 5q + 5T c p/m corresponds to the 
shift (in configuration space) between the collision points of the reference and of the 
satellite trajectories. Furthermore, 

5g c = — + - 5q c x & + a(n x <%1) 
m 2 L 

The three-dimensional maps 01]) and fl5]) remain valid for the case of planar rough 
disks. In this case all position and velocity vectors are placed in the xy-plane, and all 
angular velocity vectors are perpendicular to this plane with a single non-vanishing z- 
component (denoted Ui for disk i) [18J. Analogous definitions apply to the perturbation 
vectors. 

For smooth hard disks an elegant matrix formulation has been introduced in Refs. 
|30[ |3T] for the collision maps in the phase and tangent space. It has been generalized 
by us to the rough disk case. For details we refer to Ref. [52] . 

The equations in this section are used in the following to construct the time 
evolution of the perturbation vectors and the corresponding time-dependent local 
Lyapunov exponents. The global exponents are obtained by a time average along a 
long trajectory. 



4. Simulation results: Covariant versus Gram-Schmidt vectors 

As usual, we consider reduced units for which the particle diameter a, the particle mass 
m and the kinetic energy per particle, K/N, are unity. Here, K is the total energy, 
which is purely kinetic. Lyapunov exponents are given in units of ^JK/N ma 2 . We 
consider two-dimensional fluid-like systems consisting of iV = 88 hard disks at a density 
p = 0.7. The box with periodic boundary conditions has an aspect ratio A = 2/11, the 
same value as in our previous work on smooth hard disks [9], to facilitate comparison. 

For the computation of the covariant vectors according to the algorithm of Ginelli 
et al. [7J El [16] , the system was first relaxed forward in time for t r = 10 5 r time units. 
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Then, the Gram-Schmidt vectors were stored for t s = 5 ■ 10 4 r time units. Finally, the 
system was iterated backward in time, and covariant vectors were computed for times 
in the interval [t r ,t r + 0.9 t s ]. Here, r = 0.6 is the time interval between two successive 
Gram-Schmidt re-orthonormalization steps. 

For k > (three degrees of freedom per particle), the energy per particle is taken 
1.5 in our reduced units. For k = (two degrees of freedom per particle), the energy 
per particle is set to unity. With this choice, the kinetic temperatures agree and a 
comparison between the SHDS and RHDS is meaningful. 

The numerical accuracy for the global Lyapunov exponents quoted below for this 
system is of the order of 0.01. This estimate is based on the fluctuations of the GS- 
exponents, if the storage period r s is subdivided into ten independent time intervals. The 
agreement between the covariant and GS-exponents is of the same order of magnitude. 
Great care was taken for the choice of the simulation parameters r r and t s to allow 
the perturbation vectors enough time to rotate into position. The local exponents are 
found to be insensitive to a doubling of these parameters. The symmetries discussed 
below for the local exponents, GS or covariant, are typically obeyed with a numerical 
accuracy of the order of 0.001. For systems with only three or four disks in Section El 
these accuracy estimates are even better by at least a factor of two. This rather high 
accuracy and fast convergence is a virtue of the strongly mixing hard-particle systems 
studied here and could not be achieved with comparable numerical effort for particles 
interacting with a smooth potential. 

The left panel of Fig. CD shows the Lyapunov spectra computed with the covariant 
vectors of the RHDS for five values of the coupling parameter k e {0,0.1,0.2,0.4,1}. 
The right panel of Fig. [T] shows an enlargement for the smallest exponents of these 
spectra. Conjugate exponent pairs are plotted with the same index i on the abscissa, 
where now i e {1, 5N/2}. The open symbols are computed from the GS vectors in 
the forward direction of time, the full symbols from the covariant vectors during the 
time- reversed iteration. Considering the size of the system (N = 88), the agreement 
is excellent. The figure clearly displays the conjugate pairing symmetry A» = X^N+i-i 
expected for time reversible energy-conserving systems. 

We can see in the right panel of Fig. [1] the influence exerted by the moment of 
inertia / on the spectra. For 1 = corresponding to k = 0, the steps in the spectrum 
due to degenerate exponents are a clear indication for the presence of Lyapunov modes. 
According to the classification in [H] , the steps with a twofold degeneracy are transverse 
(T) modes - T(l, 0) and T(2, 0) from right to left. Similarly, the steps with a fourfold 
degeneracy of the exponents are longitudinal-momentum (LP) modes - LP(1,0) and 
LP(2, 0) again from right to left. Since our simulation cell is rather narrow, only modes 
with wave vectors k parallel to the x-axis of the (periodic) cell appear, leaving for the 
second argument. The modes for covariant vectors of smooth hard disks was studied 
in [9]. As usual, "transverse" and "longitudinal" refer to the spatial polarization with 
respect to the wave vector k of the wave-like pattern. As is evident from Fig. [IJ the mode 
structure in the spectra gradually disappears if k is increased. This is a consequence of 
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Figure 1. (Color online) Left panel: Lyapunov spectra for an equilibrium system 
consisting of a gas of 88 rough hard disks for various values of the coupling parameter 
(k G {0,0.1,0.2,0.4,1}). The spectra are calculated with the help of the covariant 
vectors. Reduced indices 2i/5N are used on the abscissa. Although the spectrum is 
defined only for integer i, lines are drawn for clarity. Right panel: Enlargement of the 
Lyapunov spectra for the small exponents (in absolute value). Only the spectra for 
k G {0, 0.1, 1} are shown. The open symbols indicate exponents computed from the 
GS vectors, the full symbols exponents obtained from the covariant vectors. Steps due 
to Lyapunov modes can be observed for k < 0.1, but the mode strutcure disappears 
without trace for larger n. 



rotation-translation coupling. 

4-1. Angles between the position and momentum parts of perturbation vectors 

Next, we show how differently the GS and covariant vectors span the tangent space. 
We consider, for each Lyapunov index i, the angle Oj between the spatial part and 
the translational-momentum part of the perturbation vector 5T^ - either a Gram- 
Schmidt vector g { or a covariant vector Uj. Defining the 2iV-dimensional position vectors 
$Qi — (5<li ■> ■ ■ ■ > Sq^) an d the 2A r -dimensional momentum vectors 8p i = (Spl, . . . , Sp t N ), 
we compute cosGj = • 5pj)/(||5qrJ Note, that for the reduced units 

we use throughout, the positions and momenta are dimensionless. The elements of 
the perturbation vectors consisting of the N angular velocities of the disks are not 
considered. 

In Fig. [21 the time average of cos is shown as a function of the Lyapunov index i 
for three values of the coupling parameter, k G {0, 0.1, 1}. As a reference, the respective 
result for the SHDS is also shown in the bottom-right panel of Fig. [2j Excluding the 
rotation-dominated part of the spectra (0.4 < i/5N < 0.6), one observes in the covariant 
case (red curves) and for k > that the positive branch is anti-symmetry with respect 
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Figure 2. (Color online) Time-averaged value of cos(O) = (8q ■ Sp)/(\\8q\\ \\Sp\\) as a 
function of the Lyapunov index i for a system of N = 88 rough hard disks, for three 
different values of the coupling parameter k £ {0,0.1,1}, as well as for a standard 
system of N = 88 smooth hard disks. The density is p — 0.7 and the aspect ratio 
is A = 2/11. Here, 5q and Sp are the two 27V-dimensional vectors of all position 
perturbations and of all tranlational-momentum perturbations respectively, for GS 
vectors g (blue line) and covariant vectors v (red line) . The case of smooth hard disks 
is treated in the bottom-right panel. 



to the conjugate negative branch, 

(cos (e D+ i_ i (t))) = -(cos (8i(*))) V i G {1, . . . , D/2} . 

The covariant vectors evolve according to the natural tangent bundle and obey the time- 
reversal symmetry prescribed by the respective dynamical rules. That means that the 

t h 

i expanding vector in the time-forward direction becomes the conjugate contracting 
vector —v D+ i_i in the time-backward direction. The (cos G) for the GS vectors (blue 
curves) does not show this time-reversal symmetry. Whereas the deviations for k — 0.1 
outside the range 0.4 < i/5N < 0.6 are comparatively small (bottom-left panel of 
Fig. , they are pronounced for k = 1 (top- left panel of Fig. E]). This is to be expected 
due to the fact that the GS vectors do not obey time-reversal symmetry in general as 
outlined in detail in [16] . 

It is interesting to note that the simple SHDS (bottom-right panel in Fig. |2|) shows 
anti-symmetrical behavior of (cos 0) both for the covariant vectors (as expected) and 
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for the GS vectors. For the latter, this can only be a consequence of the symplectic 
symmetry of the GS vectors for Hamiltonian systems [T6l[TT] . This was explicitly stated 
already in [22], 

cos (Q 0+1 ^(t)) = - cos (9i(t)) V % e {1, • • • , D/2} , 

for any instant of time t (which is a direct consequence of Eq. (|96|) in Sec. [7]). Of course, 
also the time-averaged spectra are anti-symmetric as displayed in the bottom-right panel 
of Fig. [2j However, the rough-hard-disk model for k = 0, keeping the angular velocity 
variables included in the GS re-orthonormalization process, does not show this symmetry 
(top-right panel of Fig. [2]). This is not a numerical artifact since the GS vectors, which 
do not display the expected symmetry, are used in the construction of the covariant 
vectors, which display the symmetry. Thus, it is not sufficient to simply put k to 
for rough disks to generate the same tangent space dynamics as found for SHDS. It is 
also necessary to explicitly exclude the redundant angular velocity components of the 
perturbation vectors from the GS re-orthonormalization procedure. 

It is well known that the Lyapunov vectors associated with the large positive and 
negative Lyapunov exponents are strongly localized in physical space [331 [23J IZE] • We 
have verified this property also for the rough disks. For details we refer to Ref. [32J. 
Not too surprisingly, the localization is even stronger for the covariant vectors, whose 
directions in tangent space are solely determined by the tangent flow and are not affected 
by re-orthonormalization constraints. The localization of a covariant vector with a large 
positive Lyapunov exponent is found to agree with the localization of its conjugate twin 
with negative exponent [32]. For orthogonal Gram-Schmidt vectors this symmetry does 
not exist. 

5. Null subspace and vanishing exponents for the RHDS 

The subspace of the tangent space associated with the vanishing exponents is of special 
interest. Let us consider some experimental results first. In Fig. [3]the Lyapunov spectra 
of a rough-hard-disk system consisting of an odd number of particles, namely N = 3 
disks, are shown. The panels in the top row are for periodic boundary conditions, which 
we exclusively used so far. The panels in the bottom row, however, are for reflecting 
boundary conditions, for which the tangent vectors change according to Eq. (20) in [19] 
for a reflection at the boundary. The essential difference is that translational momentum 
conservation does not apply in this case (both in x- and ^/-directions), which results in 
a reduction of the vanishing exponents by 4 with respect to the periodic-boundary case. 
The panels on the left-hand side of Fig. |3] include the orientation angles {9 n }n=i in the 
list of independent variables in Eq. fl2]) and the respective perturbations {S0 n }n=i * n 
Eq. (jHJ). Although the {9 n }^ =1 do not affect the dynamics in the phase space at all 
(as discussed further below), the associated perturbations {59 n }^ =1 contribute to the 
evolution of the tangent bundle since they take part in the GS re-orthonormalization. At 
the collision of particles k and I these perturbation components change [19] according 
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Figure 3. (Color online) Lyapunov spectrum of the RHDS, with N = 3 
disks, for the cases of periodic and reflecting boundary conditions, where 
the orientation angles of the disks are either included or excluded from the 
phase-space and tangent-space dynamics, n = 0.5. 
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Figure 4. (Color online) The same as Fig. |3l with N — 4 disks. 



to 

691 = 80 t - f (n x 0) 5r c 1 

59£ = 59 k -f (nxg)5r c , j {) 

which augment the Eqs. (jSJ) in this case. The panels on the right-hand side of Fig. |3] are 
for the standard representation of the RHSD without the explicit consideration of the 
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disk orientations {O n }n=i ( see ^qs. (jlj) and (E])). The spectra in the analogous Fig. @] 
are for an even number of disks, N = 4. 

The case of an even number of disks in Fig. H] is well understood whether the 
orientation angles are explicitly considered or not. In all cases, two of the exponents 
vanish due to the invariance with respect to time translation and energy conservation. 
For reflecting boundaries and excluding orientations (bottom-right panel) these are the 
only symmetries leading only to two vanishing exponents. For periodic boundaries and 
orientations excluded (top- right panel), four other vanishing exponents need to be added 
due to the invariance with respect to a uniform translation of the origin in the x- and y 
directions and to a uniform shift of the translational momenta in the x- and y directions. 
Also for the reflecting boundary case with orientations explicitly included (bottom-left 
panel), six exponents vanish, but for a different reason. The four additional vanishing 
exponents are due to the invariance with respect to a rotation of any of the four disks 
by an arbitrary angle. Finally, in the periodic boundary case with explicit consideration 
of the angles (top-left panel) all symmetries mentioned above contribute to increase the 
total number of vanishing exponents to 6 + N, where N = 4 in Fig. HI 

For an odd number of disks the situation is not clear to us. We still do not have 
a satisfactory explanation for the experimental result that an odd number of vanishing 
exponents is found, if the orientations of the disks are excluded from the considerations 
(see the panels on the right-hand side of Fig. [3]). Similarly, if the orientation angles are 
included, a naive consideration of intrinsic symmetries would suggest 5 and 9 vanishing 
exponents for reflecting and periodic boundaries, respectively. Instead, the experimental 
results shown in the left-hand side panels of Fig. [3] yield 6 and 10 vanishing exponents, 
one more than expected. One additional symmetry seems to exist, which is not properly 
accounted for in this case. 

6. Transversality versus order violation in the local Lyapunov spectrum 

The dynamics of the tangent bundle is characterized by a complete set of dissociated 
sub-bundles, which are defined by the Oseledec splitting, 



at any point T in phase space. Here, E u — Vi © • • • © ^D/2-3 and E s = t> D /2+4 © ■ • • © I'd 
are the covariant unstable and stable subspaces, respectively, and is the covariant 
null subspace (the dimension D of the phase space equals 5N for RHDS and 4N for 
SHDS). The multiplicative ergodic theorem provides us with the following inequalities 
for the global Lyapunov spectrum 



where the equal sign relates to degenerate exponents, and [A ] = is sixfold degenerate 
(assuming energy conservation and momentum conservation in x- and ^-directions). 
These inequalities, which guarantee a hierarchy of the instabilities, imply also a certain 
dynamical isolation in the tangent bundle. 



TX r = E U ® N®E S 




(7) 
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The order in Eq. may be different, if the local Lyapunov exponents are averaged 
over a finite time (instead of an infinite time for the global exponents) to which we refer 
as "loss of isolation". It may be expressed in terms of the order violation probability 
for the finite-time averaged local exponents, Aj , which is the probability that Aj~ < A^ , 
whereas for the global exponents Xj ^ Aj holds. Therefore, a spectrum is not violated, 
if all probabilities P(Aj < A i ) = for all times r > r c , where r c > 0. It implies 
that the angles between the Oseledec subspaces are bounded away from zero [31], which 
indicates uniform transversality for the respective subspaces. These considerations may 
be applied to the angles between all possible covariant vector pairs. 

In Fig. we show contour plots of P( cov Aj < C0V Aj), r = 0.01, for all j < i 
(provided A^ ^ Aj), and compare results for rough hard disks (RHDS) with parameters 
A" = 80, k = 0.1 and p = 0.7, with analogous results for smooth hard disks with 
parameters A" = 100 and p = 0.7. The averaging period r = 0.01 is small enough such 
that time averaging does not affect the local exponents much, cov Aj ' m cov A j, for all 
i. 10 5 such intervals were used for the construction of each contour plot. The regime 
UU, for which i G [2, 187] and j G — 1], refers to vector pairs from the unstable 
manifold and is located in a triangle in the upper right corner of these plots. The regime 
SS, for which % G [205,400] and j G [204, % — 1], refers to vector pairs from the stable 
manifold and is located in a triangle in the upper left corner of these plots. The regime 
SU, for which % G [204,400] and j G [1, 187], refers to the case where Vi is a covariant 
vector from the stable manifold and Vj from the unstable manifold. This regime is 
located in a square in the center of the plot, which separates SS from UU. Vector pairs 
involving at least one vector from the central manifold are located in narrow stripes in 
between these domains. A fundamental difference between the rough and smooth disk 
case may be observed. 

For smooth disks (SHDS) in Fig. [5J the violation probability P vanishes for all vector 
pairs from SU, which confirms our previous results that the minimum angle between 
the stable and unstable subspaces is well bounded away from zero, at least for finite 
A" [9]. An entanglement of covariant vectors and, hence, a loss of isolation of covariant 
subspaces due to an order violation of local exponents only occurs within the stable and 
unstable manifolds. This causes the finite numbers for the violation probability P in 
UU and SS. Our result also implies that the stable and unstable manifolds have the 
same dimension at all points in phase space. 

For rough disks (RHDS), the violation probability varies rather smoothly in Fig. [5] 
and does not vanish even for vector pairs from the central manifold. There is evidence 
that stable covariant vectors may be "entangled" with unstable covariant vectors. It 
has been argued that such an entanglement may be responsible for the disappearance 
of Lyapunov modes in particle systems [31]. We suspect that one reason for this 
disappearance may be the non-orthogonality of the central manifold to the stable and 
unstable manifold for the RHDS. 

Related and complementary information may be gained by considering the scalar 
products between covariant vector pairs. In Fig. [6] we show - for the same systems as 




Figure 5. (Color online) Violation probability contour plot for covariant local 
exponents for averaging times r = 0.01. Top panel: rough-hard-disk system (k — 
0.1). Bottom panel: smooth-hard-disk system. The order violation probability 
P( C0V Aj < cov Aj' ) is shown for all Lyapunov index pairs for which the 

global exponents obey Aj ^ A, (such that j < i). 




Figure 6. (Color online) Contour plot of (\vj ■ t>i\) for all Lyapunov index pairs 
such that j < i. Top panel: rough-hard-disk system (n = 0.1). Bottom panel: 
smooth-hard-disk system. 



Orthogonal versus covariant Lyapunov vectors for rough hard disk systems 



15 



250 

200 

^ 150 
O 

^ 100 
50 


0.000 0.002 0.004 0.006 0.008 0.010 0.012 0.014 

Figure 7. Probability distribution for the minimum angle $ between the stable and 
unstable manifolds for three rough hard disk systems of size TV as indicated by the 
labels. The coupling parameter k equals 0.4. 

before - contour plots of (\vj ■ for all j < i. The same domains UU, SS and SU 
may be distinguished. 

Let us first again consider the SHDS first (bottom panel in Fig.E]). Most remarkably, 
vector pairs from SU, which correspond to the uniformly shaded square, are not 
orthogonal to each other. This confirms our result of Ref . [9] , where we found that the 
stable and unstable manifolds are transverse but not orthogonal. The central manifold, 
however, is orthogonal to both the stable and unstable manifolds. 

The situation is slightly more complicated for the RHDS (top panel in Fig. [6]). One 
observes a very smooth variation of (\Vj ■ Vi\) for all three regimes and deviations from 
orthogonality everywhere. Most notably, these deviations persist also between a vector 
from the central manifold and another vector from the stable or unstable manifolds (the 
narrow stripes between SU and UU respective SS). To check for transversality, we 
show in Fig. [7J the probability distribution for the minimum angle $ between the stable 
subspace E u and the unstable subspace E s . It is computed from the smallest principal 
angle between the two subspaces (35j EEJ [37] . If the covariant vectors belonging to E u 
and E s are arranged as the column vectors of matrices V u and V s , respectively, the 
QR decompositions V u = Q U R U and V s = Q S R S provide matrices Q u and Q s , with 
which the matrix M = {Q U } T Q S is constructed. The singular values of M are equal 
to the cosines of the principal angles, of which $ is the minimum angle. Fig. [7J shows 
that this angle is small but positive even for the largest system considered here. The 
distributions seem to be bounded away from zero and indicate transversality for the 
respective subspaces. We conclude that for finite N the rough hard-disk systems are 
also hyperbolic, although just about. 
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Figure 8. (Color online) Test of the symplectic symmetry of Eq. (|9aj) for a SHDS 
with N = 4 (left panel) and for a RHDS with N = 4 (right panel). 




7. The RHDS is not symplectic 

It is easily seen from the motion equations that the SHDS and RHDS conserve energy. 
As a consequence, the phase volume is conserved for both systems, and 

XX*(T) = W 

i=i 

at any point in phase space. However, this is not sufficient to ensure that they are both 
symplectic. Let J denote the even-dimensional skew-symmetric matrix, 

3 = 

For Hamiltonian systems it was shown by Meyer [38] that for any Gram-Schmidt vector 
g { also the vector J'g i belongs to the set of Gram-Schmidt vectors such that 

A D G + 5 i-i(r)= -A?*(T) , (9a) 
g^J?) = ±J9i(T) , (96) 

for i G {1, . . . , D/2}. In [16] we have referred to Eq. ( 19 aj) as symplectic local pairing 
symmetry. Here, we test both expressions for the SHDS and RHDS. In Fig. |8] the 
maximum and minimum local GS exponents are shown for a SHDS (left panel) and a 
RHDS (right panel). One observes that Eq. (19 a\i is obeyed for the smooth hard disks, 
(SHDS) A G5 _ M + (SHDS) A G5 (t) = 

but not for the rough hard disks, 

(RHDS) A as Ut) + (RHDS) A p (t) Q 

A direct inspection of the perturbation vectors also reveals analogous results for Eq. 
( |9~6|) (not shown). This numerical comparison demonstrates that the RHDS is not a 
Hamiltonian system as already suspected in Sec. HI 
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The non-Hamiltonian character of the RHDS may be shown also analytically by 
explicitly considering the transformation rules in the tangent bundle [39j I2U HQ] . Any 
linear transformation S is symplectic if, and only if, it satisfies the identity 

S^JS = J . (10) 

Rewriting the linearized collision map (such as Eq.([5]) for the SPHS) in matrix form, 

ST' = SST, (11) 

we examined the validity of ( flOl) for the respective smooth- and rough-hard-disk models 
[32] . After some tedious calculation (omitted here) it is shown that Eq. (TTDT) holds 
for the SHDS, but fails for the rough disks, whether the disk orientations are explicitly 
considered or eliminated. 



8. Discussion and summary 

We study smooth (SHDS) and rough hard-disk systems (RHDS) in two dimensions, 
which resemble dense planar gases in thermodynamic equilibrium. The investigation of 
the tangent-space structure for both models reveal strong qualitative differences. Our 
approach is based on the computation of the Oseledec splitting of the tangent space 
into a hierarchy of covariant subspaces spanned by the individual covariant Lyapunov 
vectors, which are associated with respective global (time-averaged) and local (time- 
dependent) Lyapunov exponents. Wherever meaningful, we also compare the properties 
of the orthonormal Gram-Schmidt (GS) vectors with those of the covariant vectors. 
Both vector sets span the tangent space at any space point. Whereas the former 
constitute an orthonormal set, the latter generally do not. However, the covariant 
vectors reflect the time-reversal symmetry of the motion equations p2] and for that 
reason are physically more relevant than the orthonormal GS vectors [JHII2], which are 
a convenient mathematical tool for the computation of the global Lyapunov exponents. 

There is a subtle difficulty with the equations of motion for the RHDS, whether 
or not the disk orientations, which do not affect the phase-space trajectory at all, are 
explicitly excluded or included for the computation of the tangent-space dynamics (Sec. 
[5]). Exclusion and inclusion requires phase-space (and tangent-space) dimensions of 
5N and 6N, respectively, where N is the number of disks. Whereas our main results 
concerning the tangent space are unaffected by this choice, the number of vanishing 
exponents (equal to the dimension of the central manifold) is strongly affected. This 
number is determined by the continuous symmetries of the system and - if the disk 
orientations are explicitly considered - by the invariance of the motion equations 
with respect to an arbitrary rotation of any disk, resulting in N additional vanishing 
exponents. For an even number of particles, the number of vanishing exponents is well 
understood, regardless whether periodic or reflecting boundary conditions are considered 
(Fig. HI). For an odd number of particles, however, for which our experimental results 
are summarized in Fig. [3j these considerations do not give the observed numbers. This 
is a vexing problem we still do not understand. 
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Lyapunov modes, which are most familiar for smooth hard particle systems [2DJ [33] , 
exist for the RHDS only in the limit of week-coupling between translational and 
rotational degrees of freedom (Sec. HI). Already for rather small coupling constants 
respective moments if inertia, the step structure in the Lyapunov spectrum, which is 
due to the degeneracy of the Lyapunov exponents associated with the modes, disappears. 

Another observation concerns the localization of the perturbation vectors in 
physical space. Whereas the localization measure for conjugate covariant exponents 
coincides, it differs for the GS vectors (not shown, see [32] ). This is a consequence of 
the time-reversal symmetry displayed by the covariant vectors. 

As for the smooth hard disks [9J, the unstable, stable and central manifolds of 
the RHDS are transverse to each other. But the minimal angle between the stable 
and unstable manifolds typically is very small for the RHDS. For the SHDS, the 
central manifold is strictly orthogonal to the unstable and stable subspaces, which 
themselves are almost - but not strictly - orthogonal to each other. For the RHDS, 
the scalar product between two arbitrary covariant vectors does not vanish on average, 
but progressively increases, the closer the indices of the two vectors are. But it does 
not reach unity and, hence, no tangencies occur even for vector pairs with neighboring 
indices, at least for finite N (see Fig. [6]). If the global exponents in the Lyapunov 
spectrum are ordered according to size, the probability for this order to be violated by 
a pair of local exponents is always positive for the RHDS. The violation probability also 
increases progressively, the closer the indices of the exponent pair become, and may 
reach values around 0.5 for neighboring exponents (see Fig. [5]). For the SHDS, however, 
this probability vanishes, if one of the vectors is chosen from the stable, and the other 
from the unstable manifold. 

Finally, we numerically demonstrate that conjugate Gram-Schmidt vectors for the 
RHDS do not exhibit symplectic symmetry, whereas conjugate GS vectors for the SHDS 
do. This shows that the symplectic structure is destroyed for the rough hard disks due to 
the coupling of translational and rotational degrees of freedom. This interesting result 
is also confirmed by direct analytical computation. 
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